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)K>MENCLATURE 


C - Center distance 
Cg - bearing damping 
Cg - shaft damping 
CP - circular pitch 
CR - loading contact ratio 
CR^ - theoretical contact ratio 
DDELT - backlash 

DF - dynamic load factor 
E - Young's modulus 
F - gear face width 
FH - hub face width 
FW - geat web thickness 
G - torsional modulus 
GP - gear tooth pair 
HSF - hub torsional stiffness factor 
J - mass moment of interla 
JG - 1/2 Mg(RBC)^ 

K - shaft stiffness 
KG - gear mesh stiffness, N/m 
KP - gear pair stiffness, N/m 
M - mass 

P - total mesh static load, normal 

P - Hertz stress 
H 

PE - profile error 
PM - profile modification 
PSITP static angular position 
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Q 

QD 

QDT 

RA 

RATIP 

RAPP 

RABOT 

RBC 

RCP 

RCCP 

RH 

RRC 

SV 

TR 
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W 

C 

s 

y 

♦ 


static GP load, norasl 
dynaalc GP load, nonal 
total Bash dynamic load, normal 
roll angle 

RA at tip of involute 
RA at pitch point 
RA at bottom of involute 
radius of base circle 
radius to contacting point 
radius of curvature 
hub fixity radius 
radius of root circle 
sliding velocity 
transmission ratio 
deflection 
Poisson's ratio 

critical damping ratio, gear OMSh 
critical daiq>ing ratio, shafts 
dynamic displacement, rad 
dynamic velocity 

2 

dynamic acceleration, rad/sec 


Superscript: 


' - instantaneous 


V 


Subscripts : 

D • driving elenent 
G - gear 

HCR - high contact ratio gearing 
i - nesh arc position 
L - load element 

NCR - nonnal contact ratio gearing 
S - shafting 

1 - Gear 1 

2 - Gear 2 



umooMcrioii 


Many advanced tedmoloty ayplicatlona have a general rcquircttent that 
the power to transmlasion weight ratl'i be increaaed. Englneerat aa a 
result of these requlrewentSt design ;;>ar systras to mxIbub load capacity. 
However, accurate deteminatlon of gear tooth loads and stresses under 
dynamic conditions Is no^ currently possible. As a result, experience or 
engineering institution becomes the controlling factors la transmission 
design. The ability to accurately calculate the dynamic loads in geared 
systems becomes essential for advanced transmission design, 

The concern with dynamic Ipada acting on gear teeth goes back at least 
to the eighteenth century. A first concentrated effort in defining dynamic 
loads was conducted by the ASME Research Connittee on Dynamic Loads on 
Gear Te^th in the late 1920's and early 1930's. These studies presented a 
dynamic load equation more popularly known as Buckingham's Equations [1] . 

Between 1940 and mid 1950's another era In analyzing the dynamic loads 
in gear teeth developed. The studies conducted during this period utilized 
more detailed information on gear teeth deflection, and In addition, 
mass-equivalent spring models with wedge, cam, or sinusoidal type excitations 
were Introduced [4, 6, 7, 9]. In general, this group of analyses could be 
considered as using an equivalent constant mesh stiffness model. 


* 


Nundiers in brackets designate references at the end of the paper. 


Starting vlth the late 1950'a, a variable gear «eah stlffneaa aotel vaa 
considered by a nmd)er of investigators [8« 10-12* 14-16* 18-213. In these 
analyses* the gear nesh stiffness vas assuoed or calculated to be of periodic 
rectangular (or nearly rectangular) fom* in Other cases it was assuwed that 
the stiffness could be of sinusoidal or trapesoidsl foiws. The ■ein 
simplifications used either singly or in ssom coidtinations In these models 
can be generalised as follows: 

a. Gear tooth errors have negligible effect or no effect on Msh 
stiffness. This isftlies that for a given load a gear with errors 
will have equal mesh stiffness as the same gear without errors. 

b. Contact assumed to occur only on the line of action. 

c. Analysis limited to contact ratios below 2.0. 

d. The contact ratio and/or mesh stiffness is not affected by 
transmitted load* premature or delayed engagement. 

e. Dynamic simulations based on uninterrupted periodic stiffness 
functions and error displacement strips. 

In view of these limitations, the above gear mesh stiffness model can be 
defined as a fixed-variable gear mesh stiffness model (FVMS). 

The gear mesh stiffness in engagement is probably the key element in 
the analysis of gear train dynamics. The gear mesh stiffness and the contact 
ratio are affected by many factors such as the transmitted loads* load 
sharing, gear tooth errors, profile modifications* gear tooth deflections* 
and position of contacting points. 
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By introducing th«i« Mpcett, th« enlMlatnd g«ar Msh atlffnota can be 
defined aa being a variable-variable awah atiffneaa (WMS) aa oppoaed to 
the FVMS Bodeling. The need for an laiproved nodeling or variable-variable 
gear mesh stiffness modeling has been recognised or initiated to some 
degree m, IS, 16, 22-25, 27, 29, 33, 34]. 

In this study a large scale digitized approach (coiq>uter block diagram 
in Figure 1) was used for eliminating the previously Indicated serious 
shortcomings of the FVMS modeling. The concept of the WMS wes expanded by 
introducing an Iterative procedure to calculate the WMS by solving the 
statically indeterminate problem of multi-pair contacts, changes in contact 
ratio, and mesh deflections. In both the static and dynamic portions of the 
analysis, ^he gear train was ZK>deled aa a rotating system rather than an 
equivalent mass-spring system excited by the error displacement strips or 
wedges. 

The primary purpose of this study was to develop an uninterrupted static 
and dynamic analysis of a spur gear train. In both the static and dynamic 
portions, the gear train was modeled as a rotating system rather than an 
equivalent mass-spring system excited by the error displacement strips or 
wedges. At this time the modeling is limited to the condition that for a 
given gear all teeth have identical spacing and profiles (with or without 
surface imperfections). The surface imperfectlons-faults trere simulated 
by Introducing various sinusoidal profile errors and surface pitting. The 
extended modeling is illustrated by a few selected situations in the 
high contact ratio (CR^2) and normal contact ratio (CR<2) operating regimes. 
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OVERALL COMPITTER PROGRAM PLOW DIAGRAM 


1*hc dcvelo|M!i) digitized analytical method was progranNied in 
FORTRAN IV. Functional ly, the computer program is divided Into the following 
three parts: 

1. A set of subroutines to perform the static analysis. This set of 

subroutines can operate as a stand-alone unit. However, this set 
is needed to operate Sets 2 and 3. 

2. A set of subroutines to perform the dynamic analysis. 

3. A set of subroutines to perform the finite element analysis of 

gear tooth stresses. (Currently not an integral part of the entire 

system program package). 

Figure 1 depicts in general terms the block diagram for the computer 
program. The main calling program reads in and prints the input irformation 
defining the gears. It then passes control first to the static analysis section 
and then sequentially to the dynamic analysis section and to the finite element 
section. 

The MAINl routine performs the bookkeeping for the static analysis portion 
of the program. This routine calls the necessary subroutines to perform all the 
calculations required for the static analysis and the writing out of the 
results in the form of tables or X-Y plots. 

The static analysis is accomplished primarily by means of three 
subroutines: MOD, SLOWM and DEFL. 

The purposes of the MOD subroutine is to generate the XY coordinate system 
and digitize the gear tooth profiles from the addendum circle to the root circle 
for each gear. The MOO subroutine permits to build up a non-standard tooth form, 
or to introduce profile modifications, profile errors, and surface pits. 
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The SLOMM routine detenines the cmitect points and the nuaber of 
contacting gear tooth pairs, load sharing, stiffhess functions, various 
positional vectors, sliding velocity vectors, transaission ratios, etc. as 
the gear tooth pairs aovo through the aesh arc. In the SLONM routine all 
Inertia forces and torques were taken to be negligible. 

The DEFL subrmitine is used in conjunction with the SLONM routine to 
determine the individual gear tooth deflections. 

The XYPLOT routine is used to cross plot as many as four dependent 
variables against a single independent variable. 

The FAST routine is the main routine for the dynamic onalysis. The 
routine consists of a number of subroutines listed below: 

The VI BS subroutine is used to determine the eigenvalues and eigenvectors, 
of the gear train and to set the length of the numerical integration run as 
well as the integration time steps. 

The RKiriTA and the MORERK subroutines are used to numerically integrate 
the system of differential equations of motion. These routines utilize a 
fourth order Runge>Kutta integration scheme. 

The STORE subroutine is used in conjunction with the XTPLOT routine to 
generate plots of the mesh stiffness function and the dynamic force variation 
versus time. The STORE subroutine features a recirculating memory provision 
and is used as a buffer between the integration routines and the XTPLOT routine. 

The STRESS routine contains the finite element and grid generating 

subroutines to perform stress analysis of a gear tooth subjected to dynamic 

loads. At this time the STRESS routine Is not an Integral part of the entire 
program. Also see P. 74. 

The principal executing subroutines arc described In a greater detail In 


subseauent sections and Anoendlx. 


STATIC K>OF.l - mmsHMINATION OF (HIAR 
TOOnt STIFFNESS AND LOAD SHARING 
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Before the variehle-varleble aeeh etlffneae (WMS) cm he detec«taed» 
the actual (umtacting profiles auet be developed. In this process, the profile 
modificati(ms and errors eust be considered. 

It Is custoeary to define the profile aodificetions and errors by aeane 
of a profile chart. In terns of an involute chart, the profile aodifications (PM) 
and errors (PE) can be expressed as 

M ■ PV(RA) (1) 

where 

M " deviation frotn the line of action 

RA • roll Mgle limited to active profile 

PV - profile variation (error or amount of modification) as a 
function of RA 

A true Involute profile ia defined by 

M - PV(RA) - 0 


The previously discussed iJOD subroutine simulates an involute chart-gear tooth profi! 
relationship shown in Figures 2 and 5. The simulated profile chart can 
accommodate the parabolic and straight modifications of the tip and 

root zones, Figure 2a. Tne profile errors were approximated by sinusoidal 
representation. By varying the number of cycles and phase Mgle einusoldal 
profile errors (Figures 2b, 2c) could describe a large number of practical 
the theoretical cases. A sioailated surface pitting damage is shown in 
Figure 2e. The defined surface faults and their respective involute charts 
are then numerically transferred to the previously digitized true involute 


Tip and Root Modifications 
(Straight line or Parabolic) 


S 



FIGURE 2 - SAMPLE SIMULATED GEAR TOOTH PROFILE CHARTS 
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profile. This it accMpliehed by •abcraeting or oMlat tho spoelfied ooM>uiita 
of mater ial perpendicularly to the true involute profile ae ahovn in Figures 3 
and 4. 

After considering several types of "curved** segjMnts and resulting 
nuBierical difficulties, straight line aegsttnts were chosen to connect the 
densely digitised points involving the **SK>dified** gear teeth. Each gear 
tooth profile was definud by one to two hundred digitized points, depending 
on tooth size. One hundred points were used to define the gear tooth heights 
up to 25 «m (1 in.). Two hundred points %»ere used to digitise the gear tooth 
profiles above 50 on. (2 in.) In. height. The intermediate tooth heights are 
proportionally digitized between one and two hundred points. 

The digitized profile points incorporating the specified profile 
modifications and errors then are transferred to the SLOWM subroutine for 
establishing the points of contact, number of contacting gear tooth pairs, 
sliding velocity vectors, and the stiffness of the individual pairs as 
well as the variable-variable mesh stiffness. 

Figures 4 and 5 will be used to illustrate the computerized method for 
determining the VVMS and other parameters. For this purpose three 
coordinate systems are used. Following Figure 5. 

U, V - Fixed global coordinate system for the pinion and gear 
tooth profiles, gear 1 and gear 2, respectively. 

The global system, (I!, V), has its origin at the pinion 
center and its V-axis corresponds to the gear centerline. 
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X, y - Local coordinate ayatva flxad at tha root of Individual 
tooth for tha pinion and gaart raapactivalp* Tha origin 
0 (X«Y) ia located at tha intaraactiwn of tha cmtarlina of tha 
tooth and tha lino tangant to tha root circle of tha tooth. 

Tha Y**axaa coincida with tha tooth cantarlinaa. The X» Y 
coordinate ayetaai ia uaad in digitiaing tha profilao Md for 
determining the appropriate daflactiona of the teeth. 

W, Z - Intermediate coordinate ayatem rotating with the pinion and 
gear respectively. The origins of the W, Z coordinate 
systms for each gear are at the respective gear centers. 

The Z-axes coincide with the tooth centerlines. 

The transformations between the coordinate syst«ns for each considered 
gear pair (k-1, n) are: 

W1 • XI; W2-X2 

Z1 - Y1 -I- RROl; Z2 - Y2 RB02 

Ul - W1 sin PSIlTP(k) + Z1 cos PSIlTP(k) 

VI - -4J1 cos PSIlTP(k) + Z1 sin PSllTP(k) 

U2 - -W2 cos (PSI2TP(k) -1.5f) + Z2 sin (PSI2TP(k) - 1.5f) 

V2 - C- CW2 sin (PSI2TP(k) - 1.51) + Z2 cos (PSI2TP(k) - 1.5»)] 

For each angular position defined by PSIITP (k) and PSI2TP(k) the 
profile coordinates (X.Y) are first transformed into an intermediate 
coordinate systen, (W,Z), and then into a global coordinate syatem, (U, V). 

In each transformation step» the first profile point (point 1) is located at the 
addendum circle, and the final point is located at the root. Each tooth 

I 

is defined by the same mmber of digitised points. 


V 
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0(U,V) 


FIGURE 5 - GEAR TOOTH COORDINATE SYSTEM 
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For exampltt using 100 digitlssd profile points there ere: 

XI (100). Y1 (100) X2 (100) ,Y2 (100) 

W1 (100), Z1 (100) W2 (100), Z2 (100) 

01 (100), VI (100)1 Gear 1 02 (100) ,V2 (100) Gear 2 

|k ■> 1, n k ■ 1, n 

The locations of the contacting gear teeth and the number of contacting 
gear tooth pairs are determined by using a three step process. First, the 
gears are preloaded by a unit load and rotated by incrementing the 
PSIITP (k) and PSI2TP (k) angles and examining the potential contact 
between the calculated (01, VI) and 02, V2) profile points for several 
gear tooth pairs. The search technique is described in Appendix 1. The 
beginning and the end of the meshing arc are established by tracking the 
gear pair 3, (GP3) through its complete meshing arc. After the limiting 
points of mesh arc are determined the mesh arc is divided into fifty 
pcsltions. Next, the gears are fully loaded for further analysis. The 
actual load sharing and deflections are calculated for fifty arc positions by 
tracking the movement of fully loaded gears through the established mesh arc. 

By tracking five tooth pairs simultaneously, it is possible to 
analyze the mesh behavior for the contact ratios up to 3.0. Tracking seven 
tooth pairs instead of five expands the capacity of the program to analyze 
gear systems with contact ratios between 3.0 and Jtc. 

1# 

The gear tooth pair deflection 6(k)^ can be expressed in the following 

form; 6(k)^ - 6j(k)^ + fi 2 (k)j + «„(k)^ 

6j(k)i ■ deflection of the tooth of gear 1 at mesh arc position i 

" deflection of the k^^ tooth of gear 2 at mesh arc position 1 

■ localized Hertz deformation at the point of contact 



For the contacting paira, the gear tooth deflections and 62 (k)i 

Incorporate a number of constituent deflections; See Cl7] and Appendix 2, 

«1 (k)j^ “ *Ml^*^^i ^Rl^’*^^! (4) 

and sliillarly for Gear 2. In equation (A), 

5 ■ gear tooth deflection due to bending 

M 

iS * gear tooth deflection due to normal force 

» gear tooth deflection due to shear force 

.S =» gear tooth deflection due to deformation of surrounding hub 

* area (rocking action) 

= gear tooth deflection due to gross torsion of the rim or hub (Appendix 2) 

The gear tooth deflections can be considered as equivalent positive profile 

errors for the pinion and gears causing premature engagement and delayed 
disengagement [24, 29]. The presence of positive manufactured profile errors 
(material addition) will increase the total equivalent positive error at the 
point of contact thus moving it farther away from the theoretical line of 
contact and causing an earlier engagement. The negative profile errors or 
material removal at the tips will reduce the equivalent positive errors. 

In the third step, the and approportioned deflections 

were returned to equations 2 and added perpendicularly to the respective 
digitized profiles in order to simulate the above gear behavior. Now, the 
iterative search and calculate process is repeated under the "loaded and 
deflected" conditions. In this step the contacting points and the mesh 
are determined under full load. These events are illustrated in Figure 5, 
where the limiting points of meshing arc occur at points A' and B* as compared 
to the theoretical true Involute mesh arc A-B under no load. As a result, 
the contact arc, and therefore the contact ratio of the gears is increased. 

In the same procedural step, the final number of pairs in contact, locations 
of contacting points, gear tooth deflections, load sharing, stiffnesses, etc., 
are computed as the load gear tooth pairs move through the mesh arc (A'-B*)« 
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If the g«OMtricaL variations in surfaces do not pemlt contact in 
I nnd 7, tlu'O tlu* non-cont#ictlnn goar lei'tth arc Hkill subjected to 

A del loct ions. For oxaapic. If OIM and (!P3 are in contact« then for GP2 

K 

*l(2)i - ( 2 )^ and 62 ( 2 )^ - ^2^2)^ (5) 

These deflections are due to torque transmission at GPl and GP3 and 
the resulting circumferential hub defonaation at GP2. If the 
6 r 2 ( 2 )i deflections are sufficiently large to overcome the geometrical gap 

(errors) between the approaching teeth profiles of gears 1 and 2 at the 
angular position i, then the contact will be established for GF2e In 
this case the final load sharing and deflections will be recalculated on 
the basis of three contacting pairs (Step 3). These calculation methods 
can handle both the involute and non-involute gear actions, high contact 


ratio gearing, etc. 

For any mesh arc position i, the calculated gear tooth pair 
stiffness KP(k)^, mesh stiffness KG^, and load sharing incorporate the 
effects due to manufactured profile errors, profile modifications, and 
deflections by means of the Iterated numerical solutions of equations 3 
through 8 . 

The individual gear tooth pair stiffness can be expressed as 

KP(k)^ -Q (k)^ / Mk)^ (6) 

rf the effective errors prevent contact, KP(k)^ • 0. 

The sum of gear tooth pair stiffnesses for all pairs in contact at 
position i represents the variable-variable mesh stiffness KGP 

f7^ 

KOP^ - Z KP(k) 
i 1 » 

The load carried by each of the pairs moving through the mesh arc in the 


static mode can be determined as 


Q(k) 


I 


KP(k)l 

KGp, 


IP) 


( 8 ) 


vdiere p is the total normal static load carried by the gears at any mesh 
position 1 in the static node 


K 

P - Z Q(k) 

1 ‘ 


( 9 ) 



The contact ratio under non-conjugata action can ba nora proparly 
dafinad as tha ratio of tha traversed arcs* For example t referring to 
Figures and b, the loaded contact ratios for an errorless gear pair can be 
approximated as 



where 

(A* - B‘) - tPSITPl (3) - 61(3)] - [PSITPK3) - 82(3)] g, (11) 

(A* - B*)is the loaded arc length from GPS first engagMent to GPS 
disengagement with 81 being the variable angle between the tooth centerline 
and the contacting point. In this sxHlellng GP2* GPS and (?4 participate In 
the mesh arc for 1<CR<2; GPl, GP2, GPS, GP4, and GPS participate In the 
mesh arc for 2<CR<S, etc. 

For the Instances when the contact points are not on the theoretical 
lines of action (non- conjugate action) we must refer to Instantaneous pressure 
angles, instantanentis lines of action and transmission ratios. The need for 
instantaneous lines of action were indicated in [IS] and [23]. Utilising 
Figure 6, the instantaneous parameters* for the contact point A' (defined by 
RCPl and RCP2, or UCF(k) and VCP(k) in the U, V coordinate system) are: 


PPD' 

* 

RBCl/cos(a^j ■*'®B1^ 

Distance to instantaneous pitch point 

(12) 

"ai 

mi 

arcsin (UCP(k)/RCPl) 



®A2 

m 

arcain (UCP(k)/RCP2) 



“bi 

m 

arctan (RCCPl/RBCl) 



o’ 

m 

“ai * “bi 

Instantaneous pressure angle 

(13) 

TR' 

m 

(C-PPD')/PPD' 

Instantaneous transmission ratio 

(U) 

TR 

m 

RFC2/RPC1 

Involute (theoretical) transmission 
ratio 

(15) 

RBC2' 

m 

RBClxTR* 

Instantaneous base circle, gear 2 

(16) 


* Designated by superscript * 


ECCPl' -rCWPl)^ - (MCI)* 


RCCP2' -f (RCP2)* - (RBC2')* 


Equival«nt Inst ant anaous radius (17) 

of curvacura» gaar 1 


dquivaloat InatantraMia radlw (18) 

of curvatura, gaar 2 


Th« saaa procadure is used for determining the instantaneoua paraaaters as 
the hhove gear pair k traverses the aeah arc and, slallarly, for other gaar 
pairs. The instantaneous transmission ratio TR* is influenced by the deformations 
in the contact zone and tooth profile errors. It is important to note that for 
no-load and no surface fmlt conditions TR* >■ TR, and similar analogy exists 
for other parameters. 

If the actual loaded contact occurs above the theoretical line of action, 
the effective base circle radius of the driven gear will be decreased. 
Consequently, the instantaneous transmission ration, TR* will be smaller than 
the theoretical transmiaalon ratio, TR. 

In this study, it Is assumed that the instantaneous transmission ratio 
is dominated by the incoming tooth pair at point A* in Figure 6 as it moves 
through one gear mesh stiffness cycle. The approximate variation/ cycling of 
TR* is Illustrated in Figure 9. The maximum variation in TR* is defined as ATR. 

The described static analysis determines the variable-variable mesh 
stiffness (KGP), transmission ratios (TR), and the contact position vectors 
(RCPl, RCP2, RCCPl', etc.) for subsequent dynamic calculations. 
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DYNAMIC MODEL 


A gear train ahovn in Pigura 7 waa uaad in dynaalc aiaulationa. Thia 
■odtl ia aaauaad to rapraaant ona of tha practical caaaa in gearing* Tha 
nodal includaa tha input and load unite; a pair of gaara; interconnecting 
ahafta; daiq>ing in ahafting* gaara and baaringa; non-involuta action 
cauaad by gear tooth daflactiona; and loaa of contact. 

Tbe dynamic nodal is based on the sane coordinates as tha static nodal* 
The instantaneous parsMters which were detamined for various nash arc 
positions in the static analysis will also be utilised in tha dynamic 
sinulatlon* 

The equations of notion for this nodal along the Instantaneous 
(non-involute) line of action can be given In the following fora: 


D D 


♦ =BD ’d * Si ’l ♦ Ss< ’d - ’l* ♦ V ’d - ’l> 


aw 


•'oi S * Ss< S ; ♦ Ss< ’i - V * 

[CCPj(RBCl*j- RBCa'ij) + KOP^dUCWj - »BC2'Tj>] UCl - 0 

Ss'S- V 

[cGP^CRBCa'V^ - RBClY^) + KGP^ (RBC2't - RBClt^)] RBC2’ - 0 


( 20 ) 

( 21 ) 


•'l \ ^ ^ S2 ^2 ^ '^2> * *S.S < - ^2> “ ’^D * 


TR* 


( 22 ) 


T (TR') 
L 


The bearing daaping on the drive and load shafts was luapad as 


effective doiping at their respective drive and load aaseas 
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Th« braektted ttraa in equattona (20) and (21) rapraaant tha dynanic 
gear Math force which is dcpendenc on tha Oynaaiic dlsplaceaenta of engaged 
gears* gear sMsh stiff neaa and dMq>lng in tW awsh. 

In equations (20) and (21), KGP^ represents the variable^variable SMSh 
stiffnesa. KCP^ is s function of gear tooth profile errors and awdificationa, 
deflections of gear teeth, load sharing, height of engagenent, and an 
angular position 1 of engagesMnt as the gear pairs smve through the awsh 
sone* The nesh stiffness cycle Is illustrated in Figures 9 and 10* The 
basic sources of excitation for a rotating pair of gears are the variable- 
variable aiesh stiffness and the changes in Che transniasion ratio caused by 
non-involute action. The input torque is assuaed to be constant while the 
output or load torque is a function of the instantaneous transaission 
ratio shown as 1^- (TR*), ana bearing losses. Also see Appendix 5. If 

contact occurs above the theoretical line of action, the effective base 
circle radius of the driven gear will be decreased by an asKHtnC equivalent to 
the percentage decrease in the transmission ratio. 

Operational situations, which may involve momentary disengagenMt of 
gears in mesh can Impost several conditions on the dynamic gear nesh forces 
In equations 20 and 21 . By defining the relative dynamic displacement CRM as 


CRM - RBCl X - RBC2' ^2 , ^^3) 

if CRM >0 

(QDT)^ - CGP^ X (RBClfj - RBC2’y2> KGP^ x CRM (24) 

if (CRM - 0) and DDELT > (CRM| (25) 

(QDT)^ - 0 

if (CRH ^ 0) and DDELT <|CRM| 

(QDT)^ - CGP^ X (RBCl Tj - RBCl'fj) ^ (CRM-DDBLT) (26) 

Also, when KGP^ - 0, «»)T)j - 0 
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Th« •quivalmt dMKpiai in gaar Msh CGF^ vm raUtad to KC^ by aa a na of a 
critical dasploi coafficiatt (c). 



Tha Indlcatod equaClont of aotloo (aquaClcma 19 - 22) waro nuMtlcally 
Integrated In the PAST routine by aeana of a 4th order fhinge-Kutta 
Integration schose described In Appendix 5. 

The InltUl dlsplaccaents 9^(0). <*«teralned by 

statically twisting the entire systea with the prescribed T^^ and T^ torques. 
Por the Initial velocities Tjj(O), Tj(0), T2<0) and Tj^(O) the anticipated steady 
state involute action velocities were selected. 

The nuBierlcal integration of the equations of awtion is 
carried out for a length of tlae equivalent for the tia« required for the 
start-up transient to decay. This tiae Is assuaed to be equal to five tlaMS 
the longest system natural period. The Integration tlsas step Is taken either 
as one tenth of the shortest systea natural period or one percent of the aesh 
stiffness period with CR ^ 2 (two percent for CR > 2), whichever Is saaller. 
Also see Appendix 4. 

As the first step, the FAST routine calculates the dynaslc force In the 
mesh defined by equations 23-25. Next, the PAST routine Interscts with the 
SLOim subroutine to determine the adjunct dynamic Information: 

a. how the dynamic load Is shared by contacting tooth pairs during 
periods when sailtlple tooth pslrs are In contact. 

b. the variation of the load magnitude along the tooth profiles of 

a contacting tooth pair as the pair moves through the contact sons. 

c. the sliding velocity, the aaxiaua herts pressure and the sliding 
veloclty-herts pressure product along the tooth profiles. 
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FIGURE 7 - GEAR TRAIN USED IN THE DYNAMIC ANALYSIS 
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In order to save ctmputet tonal tine. It tiaa aasimed that the loaded 
neshing arcs (Points A 'and B* in Figure 6) in the static and dynanic nodes 
will be of the sane length. It is believed ttalt this is also a reasonable 
assuDption because the rapidly fluctuating loads should not produce a laating 
change of the neshing arc lengths. With this asivMption the detemined dynanic 
absolute angular displacenents can be conpared/ interpolated with the 
equivalent mesh arc positions in the static nmde PSIlTP(k), PSl2TP(k) for 
selecting the associated RCPl(k), RCP2(k), RCCPl(k), RCCP2(k) and other 
vectors for further calculations. Some of this information is illustrated 
in Figures 5 and 6« Consequently, the above-listed adjunct parameters a, 
b and c were determined by utilizing the calculated dynanic mesh force (can 
be zero for certain conditions) and interpolations bettreen the dynanic and 
static mode positions. 

For example, by establishing the correspondence between the ¥2 and 
interpolated PSIlTP(k) , PSI2TP(k) positions and associated RCPl(k), RCP2(k) 
vectors it is possible to calculate the sliding velocity for the dynanic mode. 

The necessary vector relationships for determining the Instantaneous 

sliding velocities can be seen in Figure 6. In the kinematics of gearing 

t)»e tangential velocities Vj and V 2 at the point of contact are perpendicular 

to their respective contact radii with the sliding velocity perpendicular to 

the line of action. Reference [21 ]. 

* • 

velocities and the instantaneous sliding velocity SV is determined 
by solving the vector equation 28 (vector polygon in Figure 5). 

SV(k)j « ' ^2 " - RCP2'(k) ¥2 (28) 

(Equation 27 can also be written in a scalar form as Equation 28). 

SV(k)^ - j/ (V^)^ + (V2)^ - 2V^ V 2 COS ( + a ^2> 


(29) 
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Th* dynMiic load <QD(k) for • contactiat IMr tooth pair in tha dynaaic 
mash position i was astahliahsd as 


<»)(h)^ 


^^(QDT)^ 


(30) 


For tha sama position tha Harts otraas was calculatad hy using an 
aqulvalant cylindar approach « aquation 30* 


(k), -If QD(k), 
y 1 FA 


( 


I 

RCCPl'(k) 


) 


RCCP2* (k) 

I 

whara 

RCCPl'(k)* RCCP2'(k) “ equivalent instantaneous radii of curvature 
F minimuB gear tooth face width 
A - ( ^ *^1 ) + ( ) 


El 


E2 


In this study 
(DFl)^ 


the dynamic load factors were defined as 

. (QDT)^ 

— - 


and 

(DF2)^ 


QD(k)^ 

wr 


(31) 


(32) 


(33) 


DFl can be interpreted as the dynamic load factor for the mash or as 
the dynamic load factor for the gear pair, adjacent shafts and bearings* 

DF2 is the dynamic load factor for an individual gear tooth pair 
traversing the mesh arc. DF2 is of main iBq>ortance idien the strength of 
the gear teeth is of primary importance. The larger of the two dynamic load 
factors will be defined as the dynamic load factor for design, OF. 


RESULTS AND DISOISSION 


In the extended nodeling vhlch Includes the varlsble^vsriable nesh stiffness 
(WMS) method, the gear train was modeled as a rotating syst^ excited by the 
variable-variable mesh stiffness and the profile error-induced interruptions 
of the stiffness function. The non-involute action was described by the use of the 
instantaneous line of action and consequent variations in the transmission 
ratio. 

The WMS o«thod defines the gear mesh stiffness as a function of load, 
errors and position of contact. This is in contrast with the fixed-variable 
stiffness (FVMS) method where the gear-mesh stiffness was treated 
Independently of the transmitted loads and gear tooth errors. In the FVMS 
method It is generally assumed that the mesh stiffness function is the same 
for identical gears with or without errors with the contact in both cases 
occuring only on the theroretical line of action. The non-involute action of 
the gears in the FVMS method was simulated by means of the error/displacement 
strips acting along the line of action. 

Static Analysis 

The described digitized WMS method removed many of the previous 
assumptions and simplifications thus improving the determination of the gear 
mesh stiffness. The extended modeling which includes the WMS method will 
be illustrated in the static and dynamic modes by a few selected cases in 
the high contact ratio (HCR with CR22) and normal contact ratio (NCR with 
CR<2) gearing, respectively. 

Tables lA, IB, and 2 and accompanying Figure 8 show the mesh stiffness 
characteristics for error-less gears. Presented results indicate the obscure 
but important influence of equivalent hub stiffness on the overall gear mesh 
stiffness. By increasing hub torsional stiffness (higher HSF, Appendix 2) 
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the loaded contact ratio decreaaeai Mah atlffneaa lacreaaea, chantaa in 
tranamiaaion ratio decreaaei and aenaitlvity to gear tooth errora increaaea. 

The oppoaite oecura by decreaaing the hub atiffneaa. 

The tabulated reaulta indicate aubatantial changes in the contact ratio 
with increasing loads and/or gear hub flexibilities. For exuple» starting 
with a theoretical contact ratio of 2.14 for a 32 & 96 tooth gear pair the 
loaded contact ratio can be 2.47 or higher within practical load and gear 
hub flexibility ranges. In addition* soiae of the NCR gear pairs could be 
theoretically aiade to operate in the NCR regine by selecting an appropriate 
combination of the transmitted load and gear hub flexibilities. 

Profile errors and pitting can affect the mesh stiffness characteristics 
to varying degrees. A case where only one of the meshing gears has surface 
imperfections will be considered first. With torsionally flexible hubs where 
the circionferential fixity is approximately equal to the minimum shaft size 
required to transmit the applied torques (HSFs. 5) , the sinusoidal errors of 
.013nm (.0003 in.) and narrow surface pits .Sam wide (.02 in.) were absorbed 
by the mesh flexibility without affecting the errorless mesh stiffness 
characteristics. On the other hand, when the hubs were torsionally rigid 
(HSF •> 1) the mesh flexibility was not able to absorb the errors of above 
magnitudes. Unabsorbed errors cause non-contact zones resulting in significant 
changes in the mesh stiffness characteristics (Figures 8 and 9). With 

I 

increasing hub flexibility there was a gradual return to normal mesh stiffness 
characteristics, i.e. the flexibilities in the mesh were able to narrow or 
bridge the non-contact zones. For example, a 32 & 96 tooth gear pair mesh with 
HSF « .6 was able to absorb a portion of the sinusoidal error by eliminating 
about fifty percent of the TOsh stiffness interruptions shown in Figure 9. 
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TABLE lA 


EFFECTS OF GEAR HUB FLEXIBILITY (Hi MESH STIFFNESS, 
TRANSMISSION RATIO AND CONTACT RATIO 


Gears: 32 & 96T, 8DP, 14.5**PA, F « 25.4 bw (1 ln.)» CR« - 2.U 

Nomal Load: 4450N(1000 lb ) or 175 N/nn (1000 Ib/in) 


RHl^ 

RH2. 

KG KGv 

HSF 

ATR * 

CR 

on 

f 

mm 

N^^ N/nsi^ 


X 


10.0 

14.5 

3.07x10® 1.21x10^ 

.476 

2.4 

2.47 

12.7 

18.3 

3.80x10® 1.50x10^ 

.591 

1.9 

2.42 

12.7 

38.1 

5.08x10® 2.00x10^ 

.794 

1.6 

2.36 

38.1 

114.3 

6.36x10® 2.50x10^ 

.992 

1.0 

2.32 

47.2 

148.8 

6.45x10® 2.54x10^ 

1.0 

1.0 

2.32 



TABLE IB 





LOAD EFFECTS ON MESH STIFFNESS, 




TRANSMISSION RATIO AND CCWTACT RATIO 



Gears : 

32 & 96T, 8DP, 

14.5®PA, F - 25.41MB (1 

■ in. ) , CR.|, 

- 2.14, HSF 

• .992 

Load 

»^Cmax 

KGp 2 

ATR* 

CR 


N/m 

N/m 

N/mm^ 

X 



88 

6.36x10® 

2.50x10^ 

0.8 

2.29 


175 

6.36x10® 

2.50x10^ 

1.0 

2.32 


350 

6.36x10® 

2.50x10^ 

1.0 

2.38 


525 

6. 36x10® 

2.50x10^ 

1.8 

2.41 


700 

6.36x10® 

2.50x10^ 

2.n 

2.45 



* maximum attainable stiffness in the meshing arc 

All gears without errors or modifications 

FI- F2- FHl- FH2 - 25.4nB (1.0 in.) 

RRCl - 47.2Snm RRC2 ■ 148.84nn 
^Illustrated in Figure 9 

I 


KG 


KGf- 


max 
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Gear 

Coablnatlon 

A 

B 

C 

D 


TABLE 2 

Gear Hath Stlffneaa Characteriatica 


Normal Lpad: 4450H(1000 lb.) or 175 N/mm (1000 Ib/ln) 


RHl^ 

RH2^ 

KG 

N/r* 


HSF 

ATR 

X 

CR 

10.0 

14.5 

3.07x10® 

1. 21x10 1 
2.54x10 

.476 

2.4 

2.47 

47.2 

148.8 

6.45x108 

1.0 

1.0 

2.32 

10.0 

14.5 

2.89x10® 

1.14x10 f 
1.59x10 

.714 

2.6 

1.67 

55.1 

55.1 

4.05x10® 

1.0 

1.3 

1.63 

10.0 

14.5 

3.52x10® 

1. 39x10 J 

.852 

3.1 

2.03 

37.5 

37.5 

4.13x10® 

1.62x10 

1.0 

2.6 

1.99 

10.0 

14.5 

1.63x10® 

. 64x10 1 

.517 

.6 

1.99 

37.5 

37.5 

4.89x10® 

1.92x10^ 

1.0 

0.5 

1.87 


A - 32 

& 

96T, 

8DP, 

14.5°PA, CRj 

- 2.14 

B - 20 

& 

20T, 

4DP, 

20®PA, CILj. • 

1.56 

C - 26 

& 

26T, 

8DP, 

14.5®PA, CRj. 

- 1.89 

D - 40 

& 

40T, 

4DP, 

20°PA, CRf - 

1.71 


All gears without errors or modifications 
Fl» F2- FH1» FH2 - 25.4mm (1.0 in.) 




GEAR MESH STIFFNESS, KC(N/a) 


r.PI, i:P2, CP3 


r.P2, CPI, CP4 


CP2 and CP3 


CPI exits 


CP2 axles 


IEP-1 


CPI, GP3 


GP3 enters 


CP4 enters ■ 


6P2, CP4 

Gears! 32-96T8DP, 14.5°PA 
F - 25.4SB (1 In.) 

HSF - 1.0 
Pit width ■ «5«i 
Pit depth ■ .25 » 
on Gear I* Figurc.2e 
Load«A4S0 (1000 lb) 


lost contact zones 



I 


Fj, PSITPl 


FIGURE 9 - EFFECT OF SURFACE PITS ON GEAR MESH 

STIFFNESS, (HCR GEARING) 


CEiUt MESH STIFFNESS, KC(N/«) 
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The WMS nechod can also be used to investigate other error <»^inations 
acting on both gears. For exnple, errors shoim in Figure 2b with FBI and 
”E2 of .013an are nearly setf-*coaq>enfuiting in tersis of very s«all chmgss from 
the noraal Mshing stiffness function. Other profile error coB^inatlons* 
especially of large error nagnitudes. could lead to non-operational contact 
ratios below 1.0 or to very frequent interruptions of the aesh stiffMSS 
function. The sinusoidal profile errors of approxiauitely one cycle (Figures la, 
lb) and ;01 %ib in magnitude are probably the maximum tolerable profile errors 
in accurate spur gearing applications. 

The gear tooth contacts due to deflections do not occur on the theoretical 
line of action. This results in non-involute action producing variations in the 
transmission ratio, &TR. The ATR can be viewed as a variation in the output 
torque. These variations are cyclic as illustrated in Figures 9 end 10 
could reach ST for high load and hub flexibility ranges. Some additional 
discussion of fiTR is given in section on Dynamic Analysis. 

The calculated results also indicated that the load distribution in a 

gear pair without errors remained practically the same for the considered hub 
flexibility ranges. Another observation could be made that for the gears 
with rigid hubs the attainable maximum gear mesh stiffness value remained 
approximately constant over a wide range of load. 

It is important to indicate that the FVMS and similar methods can not 
directly consider the absorbtion of errors. 

Dynamic Analysis 

The dynamic loads are influenced by a large nuober of variables such as 
the mass moments of inertia of all elements, shaft stiffnesses, transmitted 
loads, gear mesh stiffness characteristics, damping In the system, amount of 
backlash and speed. 

The presented information on dynamic loads in Figures 11, 12, and 13 ia 
intended to show the limiting ranges and effects of some of the parameters. 
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Figure 11 shove the dynsaie character let ics of a gear drive (Figure 7 ) 
vith an errorless 32 and 96 tooth gear pair, "soft" and "stiff" shafting and 
gear hubs, and varying mounts of dasping. The trends indicate that various 


gear drive system could be designed for best perfomance in tens of 
acceptable dynamic load factors DF (equations 31 and 32) by proper selection 
of masses, gear mesh and shafting flexibilities, and daiiq)ing. 

The shaft stiffnesses and the masses of the drive and load elements In 
most cases will determine Che lower natural frequencies of the system. The 
gear masses and mesh stiffness will dominate the highest natural frequency. 

The harmonic content of the mesh stiffness characteristics will excite at 
various speeds a number of natural nodes. The mesh stiffness functions shown 
in Figures 8, 9 , and in suggest a considerable variation of the harmonic 
contents for various situations. Changes in the transmission ratio TR* also 
refer to the sanse mesh stiffness cycle. The analyses tend to suggest that the 
main sources of excitation are the variable-variable mesh stiffness and its 
Interruptions. The ATR quantity which represents variations of load torque 
due to non- involute action appear to be of secondary importance as a source 
of excitation as shown in Figure lU. 

Two severe types of interruptions of the HCR mesh stiffness function 
resulting in a partial loss of mesh stiffness are shown in Figures and 10. 
The effects of unahsorbed profile surface imperfections (sinusoidal and 
pitting) are illustrated in Figures 1 '’ and 13 for the HCR and NCR gearing, 
respectively. In the presented cases, momentary gear separation can occur 
when DF>2. The resonant peaks are the average dynamic load factors based 
on the backlash between zero and .2Smm. 

The unahsorbed errors in the NCR situations considered caused a momentary 
loss of mesh stiffness resulting in hl^ dynamic loads and gear separation 
over wide regions of considered speeds. In the slow speed range 
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DYNAMIC FACTORS FOR A CHARACTERISTIC HCR 
GEAR PAIR WITHOUT ERRORS 
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there i,s a large sone of high dynamic load factors affected by a nu^er of 
the mesh stiffness function harsK>nics and separation of gears. These high 
dynamic loads can be reduced by introducing hi^er dasqping» higher applied 
loads and lower HSF*s. Reference [33] indicated a 300 percent increase in 
dynamic amplitudes caused by a aero stiffness none due to a single tooth pit. 

In this study surface imperfections were assigned to all teeth for a. given gear. 

There is a requirement for a minimum asnunt of damping to prevent the 
Mathieu-Hill type instabilities [25, 32]. In the considered cases only for 
a 32 & 96 errorless gear pair, "soft" shaft, "soft" hub and zero backlash 
case there was a narrow instability band at approximately 11,000 rpm with 

“ .05 and K “ .005. This instability was eliminated by increasing C to 

s 

about .07. The above instability could also possibly be prevented without 
changing K by including the bearing damping. However, for limiting the 
number of variables the bearing damping in this study was taken to be zero. 

There are also additional remedies for removing or minimizing these 
instabilities [25, 32]. 

The extended modeling has the capability for analyzing the distribution 
of the dynamic loads, dynamic factors, load sharing, contact Hertz stress 
(Pj^) and the contact stress-sliding velocity product (PV) for the entire meshing 
zone. The maximum dynamic loads, dynamic factors, maximum Pjj, and maximum PV 
do not necessarily occur at the same or any fixed position. These quantities and 
their locations are dependent dn the transmitted loads, speed, amount of 
damping, mesh stiffness function interruptions due to errors, and location of 
the contact points (contact point vectors and radii of curvature). Figures lU 
and 1*) show the range of the maximum P^j and maximum PV values corresponding to 
the dynamic conditions illustrated in Figures 11,1? , and 13. In general, 
these values were lower with higher damping and higher contact '•atios. 








SPEED OF DEIVUI6 







Tlie contact pregsures in the approach arc using the method of instantaneous 
radius of curvature (eqs. |H , |.<) and {l) werc« in many cases, somewhat lower 
in comparison with the true Involute solutions. The incorporated gear tooth 
deflections increased the length of the contacting and curvature vectors thus 
causing a decrease in the contact pressures. These findings are supported by 
[27 and 36]. The instantaneous sliding velocities (aq .26 ) on the other hand 
are higher than those in the true involute case. 


SUMMARY 

A large scale digitized extended gear modeling including the 
variable- variable mesh stiffness (VVMS) method was developed to analyze 
spur gearing in one uninterrupted sequence for both static and dynamic 
conditions. This approach can be used to eliminate many deficiencies of the 
currently used fixed-variable mesh stiffness (FVMS) modeling. 

In the extended modeling an iterative procedure was used to calculate 
the VVMS by solving the statically indeterminate problem of multi-pair 
contacts, changes in contact ratio, and mesh deflections. The developed 
method can be used to analyze both the normal and high contact ratio gearing 
with a minimum number of simplifications. 

The associated computer program package calculates the V\MS, the static 
and dynamic loads, and variations in transmission ratios, sliding velocities and 
the maximum contact pressures acting on the gear teeth as they move through 
the contact zone. The following findings were obtained for some typical 
single stage spur gear systems: 
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The gears and the adjacent drive and load systems con be Matched 
for optlMiM |)crformance in terns of ninlMUB allowable dynanic 
loads for a wide range of operating speeds. 

Torsional ly flexible design of gear bodies/hubs/rims offers an 
excellent means for absorbing or nininizing the ge<metrical errors 
in mesh . 

The gear mesh stiffness and its distribution are significantly 
affected by the transmit tc»d loads and tooth profile imperfections. 

The dynamic factors can be decreased by increasing the damping and/or 
contact ratio. Local damping appears to be the most efficient means 
for decreasing the dynamic load factors. 

The high contact ratio (HCR) gearing haa lower dynamic loads and 
peak Hertz stresses than the normal contact ratio (NCR) gearing. 
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APPENDIX I 


CCMTACT POINT SEARCH HETHW 


Since the location of the contact point is not constant and dlractly 
affects the amount of tooth deflection and vice-versa, it was necessary 
to develop a search technique that was able to deteralne accurately where 
contact occurred. More importantly, because of gear errors, it was 
necessary to be able to predict where contact would not happen for a 
given angular position. 

Figures A. 1-1 and A. 1-2 illustrate the search technique for the 
contacting points. In the search stage, the digitized points include the 
profile errors, modifications and appropriate deformations. Each gear 
tooth profile was described by one to tow hundred digitized points. In 
the majority of practical cases, this would trarslate into .13 to ,25 mm 
(.005 to ,01 in.) radial Intervals between two adjacent digitized points. 

Although basically the same, there are three distinct serach procedures 
in the SLOWM subroutine. These establish: 

1. Location of first contact in the meshing arc and its angular position 
(Position 1); 

2. Location of the contact point on the tooth profile as the tooth 
traverses through the contact zone (Positions 2 through 49), 

3. Location of the final contact in the meshing arc and its angular 
position (Position 50) . 

Each procedure makes at least two checks in the distance between the 
gear teeth to establish whether or not contact occurs at that particular 
angular gear tooth position. 


Tor cxampli', fur (’HlnbliHhinp, ihc ucCuAt toadud Initinl point ol 
the meshing arc. Point A', Figure 6, the geers vere counter rotated to 
be outside the theoretical initial contact point equation A. 1-1> Then 
the loaded gears are rotated by subtracting snail Incrments (DELTA x 
NLIM X DPELT) from PISL in the direction of actual rotation and conparing 
Che gaps between the approaching loaded teeth. 


PISL • PSSISL + DELTA x (NLIM x DELT) 


(A. 1-1) 


PISL 

PSSISL 

DELTA 

NLIM 

DELT 


starting search angle 

theoretical starting angle for meshing arc 
angular Increment In the contact aone 
arbitrary number such that NLIM x DDELT ■ any integer 
greater than 5 

1/N, where N Is any interger greater than cne 


The product (NLIM x DELT) determines how many angular increments, 
DELTA, the gear teeth are set back. The product (DELTA x DELT) 

DELTA, the gear teeth are set back. The product (DELTA x DELT) will 

later determine how much the gear teeth are incremented in the search 

for the Initial contact. Generally, for light loaded, fairly rigid, 

non-modlfied gear teeth, (NLIM x DELT) can be In the range of 5 to 8. 

NLIM should be larger for systems with modified teeth or relatively 

soft teeth or hubs. Smaller values of DELT will give more accurate 

results as to the location of the initial contact. Both large values 

of NLIM and small values of DDELT result in more InCerations, more 

accurate results, but require more computer time. 

Each gear tooth is described in space by the U1 (J) , VI (J) and 
(L), V2 (L) coordinates. In searching for the initial point of 
contact, the search is stirted with the tip of the driven gear, point 
L ■ 1 by examining the gaps between the tooth profiles of both gears. 
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Figure A. 1-1 


SEARCH FOR INITIAL MESH ARC CONTACT 
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For any 

VI (J) i V2 (L) > V1(M) 

nil can bo calculiitod by Htnllar triangles fro* Figure A. 1-1 to be 

If. Ull > 02 (L) 

and 

(Ull - 02 (D) S 0.00010 in. or .0025 on., (A. 1-3) 

a penissible aaount of Janaing or overlap has occurred; contact is 
established and the rotational angles for this position defined (Point A' » 
Figure 6). 

If the condition (A. 1-3) was not satisfied by any of the digitized 
profile points, then the sane points were reanalyzed for the "ainimum 
contact/gap" condition (A. 1-4) 

absjuil - 02 (L)| (.00001 in. or .00025 nn.) (A. 1-4) 

Using this second test, a contact was declared if the condition 
(A. 1-^) was satisfied. 

There is no contact, if neither of the above two conditions are met. 
In this case, L is incremented by one, and the process repeated. If all 
the values of L in the search region have been exhausted and no contact 
found, then the angular position of the gears is advanced by an amount 
(DELTA X DELT), and the search process repeated with L ■ 1, etc. 


The technique for finding the final contact point is similar to 
the one Just described. This time, the search is initiated with the 
tip of the pinion point J, J • 1. 

After the initial and final contact positions are found, the 
contact positions for the remaining 48 meshing arc positions are determined. 

For any i^^ mesh arc angular position, the contact points or absence 
of them are established by analyzing from 20 to 40 digitized points for 
each gear pair profiles in the approach or estimated contact zones. This 
is acciMsplished by incrementing the vertical search distance Vll (common 
for each pair profiles) and comparing the corresponding horizontal "U” 
distances between the profiles by means of the previously discussed 
conditions A, 1-3 and A. 1-4. 

The allowance (A, 1-4) was introduced to account for small deviations 
in the profile digitizing and other numerical processes. It should be 
noted that the longer horizontal rather than the shorter perpendicular 
distances were analyzed thus increasing the probability of contact. The 
(A. 1-3) and (A. 1-4) conditions were established by investigating a 

large number of gears for the known theoretical "contact" and "no contact" 
points. For the situations falling both tests, there was an unacceptable 
gap or no contact. In the search method the initial or the highest point 
for the pair in the i^^ angular mesh arc position. 

This process was repeated for all tooth pairs expected to be in 
tH 

contact in the i mesh arc position. Referring to Figure A. 1-1, contact 
was established at point P for the gear pair GP (k+1) . There is no contact 
for GP (k) . Next, the gears are advanced to a new angular position, and the 
process is repeated for the entire loaded mesh arc. 
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DEFLECTIOl AT POINT OF CONTACT 


Nunerlcal lntet;ratioR of digitized gear tooth slices (Figure A. 2-1) was 
used to obtain the bending (A ), shear (£ ) and normal force (5„) deflections 
used in equations 3 and 4. These calculation were performed in the DEFL subroutine. 
The circumferential deformation of the gear hub and deformation of the 
adjacent part of the gear body were reflected to the contacting point as 
(5 ) and (4_) deflections, respectively. 

The methods for calculating the 6 and the localized hertzian deflection 

o 

6 are amply described in [3, 17, and 29]. The 6 , , 6 and 6 deflections 
H R1 B H 

were calculated in the SLOWM subroutine. 


The 6^ deflections cannot be easily defined. Following [17]* these 
K 

deflections can be approximated for Gears 1 and 2 as shown below. 




Q(k)i (RCPl(k)i) COS n B1 r . 
4 It Gl(FHl) '■ ^ 


1 / , 1 
“rhit ^ 


-) ] (A. 2-1) 


where 


Q(k) , * load along the instantaneous line of action at the 

^ th 

contact point, k pair. 


RCPl(k)* radius to the contacting point, Gear 1, k^'' pair 

FHl » hub face. Gear 1 

RHIq » outside hub/rim radius, Gear 1 

RHlj * effective radius of circumferential hub fixity. Gear 1 

G1 =» torsional modulus of elasticity. Gear 1 


Similarly, for Gear 2 

&.,(k) - ^ ^CP 2 ( .k >ij . ^ cos g B^r/ i 1 )^] 

4 - G2fFH21 RH2. ^ ' RH2 ' ^ 


(A.2-2) 


In many cases, It could be assumed that and will be 

approximately equal to RRCl and RRC2, respectively. The radius of circumferential 
fixity for individual hubs cmmot be as readily assumed. will depend 

on the hub disk face width (>1F), hub web thickness (HW), type of gear 
mounting, shaft size, cutouts, etc. 

The torsionally rigid htibs can be theoretically obtained when the 
radius of circumferential or torsional fixity will coincide with the root 
circle resulting in 5 > 0. The opposite case can be visualized with the 

thin hubs being fixed to small shafts. In general, an increase in the hub/rim 
flexibility will Increase the total deflection of the tooth and thus will 
decrease the gear mesh stiffness. The hub stiffness factor (HSF, eq. A2-3) 
will be used to indicate a degree of Influence of the hub flexibility on the 
overall gear mesh stiffness. 


HSF 


KGmax 

KG*max 


(A. 2-3) 


where 

it 

KG » maximum mesh stiffness with torsionally rigid hubs 
max ® 

KG * maximum mesh stiffness with designated hubs 
max 

A combination of rigid hubs will be identified by HSF » 1. 
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PROFILE DIGITIZING CUBROUTIME-MOD 


The MOD aubroutinc ia used for di^itizinf; the spur Rear tooth 
profiles • Both standard and non-standard gear forms can be digitized. 

This subroutine can accomodate the parabolic imd straight line modifications 
of the tip and root zones, sinusoidal profile errors and surface pits as 
shown in Figures 2 and 3 and discussed in the static analysis section. 

The main parameters needed for describing standard and/or non- 


standard profiles for each gear are: 


DP 

M 

PHID 

TG 

AD 

WD 

GRRF 

PATM 

STTM 

RATM 

PABM 

s’i:m 

RABM 

PER 

PAP 

CYC 

IPIT 

DEEP 


diametral pitch (English input only) 

gear module (metric module only) 
pressure angle, degrees 
number of gear teeth 
addendum 
working depth 

generating fillet radius of basic rack 

parabolic tip modification 

straight line tip modification 

roll angle of tip modification, degrees 

parabolic bottom modification 

straight line bottom modification 

roll angle of bottom modification, degrees 

amplitude of sinusoidal error 

phase angle of sinusoidal error 

number of cycles of sinusoidal errors 

profile coordinate points over which pit occurs 

depth of pit 


Other symbols used in the computerized profile equation in the MOD 
subroutine are defined in the Program Listing, Apnendix B. 

A number of figures are included in this Appendix to show the 
graphical relationship of the principal profile-defining symbols. 

Figure A. 3-1 shows a basic standard involute tooth profile. Figures 
A. 3-2 through A.3-L depict several modifications of a standaid involute 


tooth profile. The fillet radius RFl for gear 1 is described as 

RFl = .7*{GRRF1+WD1-AD1-GRRF1*»2)/ 

(.5*PDl»WDl-ADl*a3FFl) (A. 3-1) 
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imd, sinilarly, for Gear 2. In equation A. 3-1* PDl ia tUe pitch diawetor. 

» 

Other symbols are defined at the beginning of this Appendix. 

The sinusoidal profile error PEI for Gear 1 was defined as 

PEl(J) ■ PER1«SIH((»RATI1-RA1(J)) 

*CYCl/RATIPl)+PAPl) 

and, similary, for Gear 2. 

In equati,on A. 3-2 the phsise angle PAPl refers to the peidc of the error 
from the pitch point. The sinusoidal error covers the region between the 
tip and root nrofile modifications as shown In Figure A. 3-3. 

A straight line tip and root profile modification model in terms 
is shown in Figure A. 3-3. A similar model was used for Geai* 2 as well 
as for the parabolic tip and bottom modifications. 

By introducing negative profile modification in the root zone as shown 
in Figure A. 3-^ several types of undercuts can be developed. 

The program has several protective features. For example, in the 
case of very severe profile modifications the contact ratio could fall 
below 1 or there could be an interference, then a special notice will 
be printed and the program execution stopped. 
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Case: RBC > RTF 

RBC > RRC 
RLM > RTF 



RRO 

I 


! . 



RTI, RAC 


BETA 


nj - Normal to profile 
6j - THETA (J), Normal angle 
RTI- Radius to top of involute 
RAC- Radius to addendum circle 
RPC- Radius to pitch circle 
RBC- Radius to base circle 
RLM- Radius to limit circle 
RRC- Radius to root circle 
RLI- Involute limiting radius 
RF - Filled radius 
RTF- Radius to top of fillet 


^O(W.Z) 


Figure A. 3-1 


PRINCIPAL SYMBOLS - 

STANDARD INVOLUTE TOOTH PROFILE 
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N. normal to original 
^ involute 

N normal to modified 
° profile 


THETA(J) 

THEM(J) 


■•■PE material addition 
-PE material deletion 


Figure A. 3-2 MATERIAL ADDITION OB SUBTRACTION FOR 
A TYPICAL PROFILE LOCATION 


*1 



RAl(J) 
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Figure A. 3-3 STRAIGHT LINE MODIFICATION OF PROFILE AT TIP AND ROOT 
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end of involute « RAH 


end of original involute, RAB 


Figure A , i-U 


TOOTH ROOT MODIFICATION-MATERIAL REMOVAL 



APPENDIX 4 


VIBS SlIBHOUTINE 


The VIBS subroutine based on a Jacobi interation technique is used to 
determine the eigen-values of the gear train from which the length of the 
numerical integration as well as the integration time step were determined. 
Using the undamped version of the model shown in Figure 7» the equations of 
motion expressed in matrix form become 

[J] {?} t [K] {V} - {0} (A. 4-1) 


The inertia matrix (J) is 

|- 

*DJD 0 0 0 

0 DJGl 0 0 

0 0 DJG2 0 

000 djl 


(A. 4-2) 


The stiffness matrix (K) is 

dktif: -dkds ° ° 

-DKDS DKDS + DKAVG x (DRBCl)^ -(DKAVG X DRBCl x DRBC2) 0 

/ 2 
0 DKAVG X DRBCl x DRBC2) DKDS DKAVG x (DRBC2) -DKLS 

0 0 -DKLS DKLS 

— ^ 

l)K'AV(* is till* mcsli stiffness. It Is determined by summing 

up the stiffness function over one cycle and dividing by the number of 
positions in the cycle. 


1 

DKAVG - - r KCP. ; n - lEP - 1 

" 1-0 ‘ 

* 


(A. 4-4) 


Shown in double precision format 
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In equation A. A-4: 

KGPj^ - gear i^fih stlffneae at the 1^** position of the nesh arc. 

n number of the neah arc positions in one stiffness cycle 
lEP ■ start position index for new stiffness cycle (Figure 9) 


Having defined the constituent parts of the above utrix equation, the VirS 
routine is called to determine the natural frequencies and the oiodal shapes. 
The natural frequencies were used to determine the tine period o.er which the 
system is to be evaluated and the length of the integration tlt^^a steps. It 
was assumed that the startup transients ifould decay within a tisM period 
equivalent to five times the longest natural period. This time period plus 
the time required for 2 or 3 additional cycles, depending on the contact 
ratio, constitutes the total time span TTOTAL in the Integration portion of 
the FAST routine. The integration time step DT is taken either as one tenth 
of the shortest system natural period or one percent of the stiffness period 
with CR <2 (t«m percent for CR > 2), whichever is smaller. 



APPENDIX 5 
PROGRAM INTEGRATION 


Thr routine hur. been developed to analyze a four macint mathematical 

model of a neared torsional system shown in Figure 7> This model 
also includes the gear mesh stiffness variations i damping in the shafts* 
gears and bearings; non-involute action effects and loss of gear 
tooth contact due to dynamic conditions. 

The Fast routine (Figure l), which constitutes the dynamic portion 
of the entire progr/un, consists of the natural frequency subroutine 
VIES, the integrating subroutines RKUTTA and MORERK, auid the storing 
and plotting subroutines STORE and XTPLOT, respectively. 

The differential equations of motion were programmed in the double 
precision on IBM 370/158 computer. The eqtiations of motion (equations 
19 - 22 ) are shown in the ccanputerized symbols as equations A. 5-1 
through A. 5-^. All variables used in these equations were declared 
double precisioji. Thin also includes the symbols preceded with letter D. 


PSDPUU = (-dcds*((psdpd+d(m;ad)-(psipd+domgai)) 

-DCbjJ*’? ( PSDPD+DOMGAI^DCBl* ( PS1PD*D0MGA1 ) ( A . 5-1 ) 

-DKDS* ( ( PS0P+DT*D0MGAD ) - ( PSIP+DT^DWGAI ) ) + TD ) /PJD 


PSIPDD = ( -nCDf • ( ( PSIPD+DOMGAI )-(PSDPD+D0MGAD ) ) 

-DKDS» ( ( PS1P+DT*D0MGA1 ) - ( PSDP+DT*DC»<GAD ) ) ( A . 5-2 ) 

-CGP« ( DRBC 1» ( PSl PD+DOMGAl ) - ( KRBCN»PS2PD + DRBC2*D0MGA2 ) ) 
•DRBCl 

-KGP*(DRFIC1»(PS1P+DT*D0MUA1)-(DRBCN»PS2P DRBC?»OT 
«1/0>«GA2) )»DRBC1) /DJGl 


PSPPDD » (-DCLS«((PS2PD+D(»«JA2)-(P.SLPD+DOMGAL)) 

-DKLS»( (PS2P+DT*DaiGA2)-{PSLS+DT»DOMGAL) ) (A . 5-3) 

-C0P*( (DRBCN»PS2PD+DRBC2‘*D0MGAn)-DRBCl»(PSlPD+D0MCAl) ) 

•DRBCN 

-KGP« ( ( DRBCN»PS2P»DT*D0MGA2 ) -DRBC1» ( PS1P+DT»D0MGA1 ) ) 
•DRBCN/DJG2 

PSLPDD « (-DCLS»((PSLPD+D0MGAL)-(PS2PD+DCMGA2)) 

-DCBL*(PSLPD+D0MGAD- DCB2*(P52PD+P0MGA2) (A.5-1i) 

-DKLS* ( ( PSLP+DT«D0MGAL ) - ( PS2P+DT»D0MGA2 ) ) - TL ) /DJL 


$5 


effeetive le*d torque It equetion (A. 5-^} it teorouwA to teeotatt for tbo 
botrlng lotMt Iqr redefining the loed torque nt 

TL ■ (TI>-DCBD»(PSDFIH-DONSAO) « DCll*(PSIFIHl)(liaAl})«1iai 
- 1)CB2»(P82FIHD0NQA2} - DCM.*(PAIfIHD0li3AL} 

For maintaining greater numerical accuracy Iqr working with laz^r numbera the absolute 
angular displacements and velocities were introduced into equatlcMis through 

A.5-^« For example, in equaticm A.5-2t DPSIP is the angular oscillatory displacement 
and DT*DQMSA1 is the swept out constant angular displacmsent of Gear 1. ^e absolute 
angular velocity at any instance consists of DPSIFIH-IXMQAI terms where DPSIPD is the 
oscillatory component of the constant angular velocity of gear 1, DOMGAl. Similar 
expression were introduced for gear ?. by using the effective base circle radius 
DRBCH, DT and C0M6A2 values. Ihe initial displacements were determined by static8Q.ly 
twisting the entire system with the drive and load torques, TD and TL, respectively. 
Thus, at time equal zero, the initial displacmaents are: 


Text 

F^(0) 


C(xnputer ProKram 

DPSID » TD/DKKS (radians) 


- 

DPSn 

= 0.0 

'I'gCo) 

- 

DPSI2 

a - TD/KGP X DRBCl x DRBC2) 


- 

DPSIL 

= DPSI2 - TL/DKLS. 


The initial velocities are set to ncmiinal steady-state velocities 


Text 


Computer Prcgram 


yo) 

’ »j(0) 

DPSIDL » 

BPSIID » 

2 x RPMIN/60 (rad/sec) 


' »i(0) 

DPSI2D a 

DPSTLD ■ 

(DRBC1/DRBC2) x DPSIDD. 
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n»e general integration scheaetic for pro<’ra« integration In the FAflT roatine ia 
shown In Figure A. 5 - 1 . 'iije fictiml numerical inte^ation la perfors»d In the 
RKUTTA and MORERK subroutines baised on the fo^irth order Bui^e~Kutta method [2€J« 

The RKUTTA subroutine (Figure A. 5-2) keeps track of the iteraticms across the 
integration interval. The MOWSOC subroutine (Figure A. 5-3) evaluates the 
derivatives and performs the suasMtions. 

The RKUTTA call statement arguswnt contains the integration step sise. The 
MORERK ccdl statement argument contains the variable to ^ integrated » its 
derivative value, and the Integration time step. At first glance, it would se«n 
that the variables are being integrated in reverse order. But, it must be 
remembered that the integrated values asre those that will be used in the next 
integration step. Therefore, the positional values are integrated first, and then the 
angular velocities. The MORERK subroutine is called eight consecutive times after 
RKUTTA to evaluate each element's change in position and velocity; PSDP, PADPD, 

PSIP, PSIPD, PS2P, PS2PD, PSLP, and PSLPD. The variable NE is a counter used to 
index the Integrated variable and its derivative in two- eight element vectors, 

XI and DXI. It is reset to zero at the start of every interation. HP is the 
vnrinble controlling the iteration time step and denoting the iteration step for 
MORERK. NRK is a V6u*iable used in the calling subroutine to check for the conclusion 
of in'-'>gration for a given time interval. 

The integrated values of angular displacements and velocities (equations A. 5-la 
through A.5-lta) for each element represent the deviations from the nominal constant 
velocities emd swept out displacements. These values are when added to the constant 
velocities and swept out displaceiMnts to give the espective absolute angular 
velocities emd displacements. In addition to being used to initialise the next 
integration step, the absolute angular position is used to interpolate a new value 
for KGPj and TRHj^. 
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Caletilntlou in FAST «r» bnsnd «» a itlffhiM eyeUae shotm in Figure* 8 »nd 9< 

!Itoe cycle starts eith the Initiation of eotttaet os a tooth entering the contact 
sons and ends with the initiatiMi of e«Btaet vith the tooth iipMdiately follovlng 
it. In the prograa, this is done in the 8LQIM si^routlne hy exaaining the developed 

^ a. 

stiffness functitm. The position of tooth #3 vben coaes into eontaot is defined 
as lEP. Consequently, (IEP-1) is the endpoint of the stiffness cycle started idien 
tooth # 3 cane into contact as illustrated in Figure 9* This process is repeated 
until the total number of stiffness cycles (IK^) equivalent to TTOTAL is reached. 

At this point, it is assumed that the systen is at a steady state and the valiMS 
of gear pair stiffness, dynamic force, angular position, angular velocity, 
stifihess, hertz stresses and dynamic load factors are printed out. The RCPl, 

RCP2, RCCPl, RCCP2 vectors needed in scane of these calc\ilations were similarly 
interpolated as KGP. 

The total integration time TTOTAL and the integration time step DT are based on 
the lowest natural frequency of the systes. More details on TTOTAL and DT are 
given in Appendix 

The description of parameters and relationship between the text and conputer 
program symbols is given in Table A, 5-1. 



Figure A. 5-1 PROOtAM IHTBGRATIOR 




Figure A. 5-2 SUBROUTINE RKUTTA 

















fmE A.5-1 


RBLASI0R8RIP BBMSil TOR 

fsa m mstam emm 


TEXT 

PROQRAN 

i»8CRzmo« 

V 

PR 

Poleson't ntto 

*1.*2 

TDEPLl, TOma 

Toeth defXeetlM, Oeerlu^ 2 

C 

ZKFAO 

OMr>aesh« erltleel dioqpiiH( ratio 

Cs 

ZITAS 

Shaft oritioal dopiiig ntio 

TR* 

TRX 

Inetantaaeous traasaiaalon ratio 

BBC2* 

»BCH 

Instantaneous bause circle, Gear 2 

RCCPl* 

RCCP2* 

RCCl 

RCC2 

Instantaneous radius of cur-rature 

Cj3 

[»•<] 

Inertia matrix 

m 

[3M] 

Stiffness matrix 

»2 ^2 '^2 

PSK>,PSDPD,PSDPDD 

PS1P,PS1PD,PS1PDD 

PS2P,PS2PD,PSPDD 

PSLP.PSLPD.PSLPDD 

Dynamic Displacement, Velocity, 
and Acceleration, Driving Gnit 
Dynamic DisplaceaMnt, Velocity, 
and Acceleration, Qe^ 1 
Dynastic Displacei^ct, Velocity, 
and Acceleration, Gear 2 
Dynamic Displacement, Velocity, 
abd Acceleration, Load Gnlt 


TD,TDIM 

TL,TOUT 

Input torque 
Output torque 


OMGAD 

Constant angular velocity, driver, 


omkx 

Constant angular velocity, gear i.rad/sec 


0MGA2 

Constant angul.jx velocity, gear 2, rad /sec 


OMGAL 

Constant angular velocity, load, rad/sec 


NOTE: Letter D jareceding the cooputerlzed tyaibols identity doxible precieioo 
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PROPOSED PUITB EUMait PROOtRN « STRESS 

IttTROlHKmOH 

The Mternel spur gcer progm lAich mm developed by the Kec^leel Bngipeering 
Department, Cleveland State Univeralty tmder NASA grant NJ^-ISSA? currently liaa 
Che capability to aolve for the inatantMooua atatic Md dynamic lo^a of an 
external spur gear mesh. The program alpr> calculatea the inatantMeoua tranamlaaion 
ratio, eliding velocities, naxiawm cMtact pressure, Md maximum FV values. For 
determining the gear tooth bending stresses the finite elMMt approach is 
proposed. 

Several types of elMients and mesh generators were considered. In order to 
conserve computing time without compromising the accuracy, it was decided to use 
an isoparametric inccmpatible displacement finite eleMnt. The entire Mternal 
spur gear program was developed in a modular form, thus facilitating introduction 
(if needed) of other types of finite elements or even the major finite element 
programs such as NASTRAN, Sx\P IV etc. 

The quadrilateral elements used are a four node quadrilateral |uiving two 
incompatible displacement modes that were developed by Wilson, et. al., [37]. 

This element was selected because it gives excellent results for bending 
applications, while maintaining a narrow band width. For large computer systems 
where the band width is no problem, the more common eight or nine node 
isoparametric elMcnts could be used. The incompatible displacMent model is 
discussed in more detail in the subsequent sections of the appendix. A special 
purpose finite element mesh generation routine was developed that takes into 
account the digitized spur gear tooth shapes from the HOD routine thus minimizing 
the input on the part of the analyst. The static and dynamic portions of the 
external spur gear progroa contain the required definitions of the gear tooth 


gMiwtry, definitions of the Mstsrial properties es estl se the ImfMiMMMNNM " 
gear load» its tocet^ «id orientation. 

The finite alsMat si^ aeah geoaration rbutteea HMra cattlaad aa a sla|la 
STIESS routine as ahoen in Figure A. 6«1. Linking* the proposed STkiSS ' 
routine with the shove described external spur gear progroi la a function of 
coi^uter alia Md type. This finite elanent prograa* due to its lar^ scale* 
at this tine la not an integral portion of the gear load analysis progran. 
Currently the proposed STRESS progran uses the pre»procesaed data from tiM 
above external spur gear progran - there is no direct Hide het«>een the load 
and finite stress analyses. Work la to be continued to develop a direct link 
betareen these t%io programs for usage on various c<Msputers. 

In the STRESS program the stresses can be calculated on a gear or pinion 
tooth at any profile position along the line of action. For efficiency* the 
stiffness eatrix is assembled and decomposed only once for each structure 
as shovn in Figure A. 6-1 . 

The decomposed matrix Is then stored on a disk file for subsequent use for 
each selected load case. Since oiost of the computer time involved in the 
stress calculation is in the assembly and decomposition procedure* this 
approach represents a considerable savings in execution time for fear sets 
having the stress routine called more than once. The following load cases 
can be called: 

a. Maximum dynamic load (at any position). 

The calling index FELGR ■ 1* Control Dame List* CONTRL. 

b. Dynamic load at pitch 

The calling index FRLHR •• 2* Control Name List* CONTRL 

c. Maximum dynmlc bending moment 

The calling index FELGR > 3* Control Name List* CONTRL 

d. Maximum static load (at any position) 

The calling index FELGR ■ 4* Control Name List* CCNTRL 


mBummm 

This pregrM takes ^vsittst*' of generic eh^ of the qnir goM tooth. , 

The finite eleeent eesh is divldi^ into two s^srste srees* the tooth (figure A.6-S) 
end the geev blank (figure A. 6~3) nodeled as a triaagulsr seme. 

The gear tooth portion is coi^eed of quadrilateral aleeenta as illustrated 

in figure A. 6-2. An algorittoi vaa devised that CMsea nodal points near the 
gear surface to be aore closely spaced than those near the centerline of the 
gear tooth. In the horizontal direction using syeMtry six nodal points 
tiere selected. The relative coarseness in the horlstmtal direction - X 
direction is shown in Figure A. 6-2. In the vertical directlmit the total 
number of nodal points isust be even. It la felt that about six to eight 
nodal points will be sufficient to describe the working portion of the profile. 

The fillet zone can be described also by 6 to 8 nodal points. The vertical 
distance between two adjacent nodal points encompass several digitized profile 

points. The relative coarseness between the nodal points on the gear tooth 
profiles is suggested in Figure A. 6-2. 

The only input that is required for the stress analysis is the desired total 
number of the surface nodes on the gear tooth profile (MNODE) and their respective 
index numbers (MODES). The name list heading is FIMLEM for HMODE and IWDES. The 
coordinates of the desired nodal points are contained in the NOD and STRESS 
routines. The grid generator automatically generates the required data from the 
designated surface nodal points for the analysis of the tooth based on a sysmetrical 
side input. Coordinates of the involute profile points are transferred by a 
common block statement. The mesh is generated by defining the masher of nodal 
points (NNODE) and their respective index number (MODES). 

The lower portion of the structure is a triangular shape as illustrated in 
Figure A. 6-3. Point 0 is the center of the gear. This area is composed of 
triangular and quadrilateral elements. The locations of the nodal points are 
dictated by the nodal points at the base of the tooth and the geimetry of the 


•f trlwslt AOB M llluttr«t«d. No «o4«l point ii^tn «r« roqwirod in 
this CMO. 

IMCOHPATXBU DlSPLACBfBIT 
FINITE BLEMBIIT 

This nnslysls uses sn iTOpsranstric finite elssMt hnvlng tno Incos^tiblo 
dlnplncsAMt nodes. This type of elesMnt ellowe for e lioMr etrein field 
heving cos^reble eceureey for bendlnp problems of eight or nine node 
ieoperenetric elements* but result in e bend width Identicel to e constent 
strein element. 

The basic probloa with constant strain elonuits is that they behave poorly 
tmder pure bending. When using this type of element for the analysie of a 
structure in a location when bending stresses predominate (for exMiple 

in the area of a gross structural discontinuity on a pressure vessel) the 
stress analyst must take special care to define a finite element mesh that is 
sufficiently refined to assure reasonable convergence. As a result* this type 
of analysis often requires large amounts of the stress analyst's time to 
be devoted to setting up analytical models and to interpret a large amount of 
output data. Also, by refining the finite element mesh, an increased number 
of equations must be solved. This may appreciably increase the computer time. 

The usual approach to circumvent these problems is to use a higher order 
element. This type of element uses a higher degree polynomial to approximate 
the displacement field and usually has one or more mid-side nodes. One 
type of higher order element idiich has no mid-side nodes and is used in 
this study tias reported by Wilson, Taylor* Doherty and Ghaboussl [37]. 


The approach uaad by VllaoBt at. al.\ In thair bi^r orter alaaant 
davalo|«ant la to add hl^r ordar tanu to tha diapl^^Mit nodaa of loaar 
ordar alaiianta to conpanaata for tha arrora In tha lowar ordar almant. A 
brief dlicuasloB of their technique follotra. 


Conaldar tha two-dlnanalonal alenent uhleh la llluatratwl In Figure A. 


The exact diaplac«Mnta are of tha font 


u ■ c. xy 


V - 1/2 Cj (a* - x^) + Cj (b^ - y*) 


(A.6-1) 

(A.6-2) 


Theae dlaplacenenta are illuatratad In Figure A. 6-4b« For a conatant 


atrain eleaent» the aaatottd dlaplacananta are of tha fom 

u - Cj xy 


V ■ 0 


and illustrated In Figure A. 6-4C. 


(A. 6-3) 
(A. 6-4) 


From equations A.6-1 through A. 6-4, it can be seen that the error in the 
constant strain solution is of the form: 


V ■ (a^ - x^) + dj (b^ - y^) 


(A. 6-3) 


Wilson et al. proposed that the error described by equation A. 6-3 could be 
eliminated by adding displacements of the form of equation A. 6-5 to equations 
A. 6-3 and A. 6-4. In this way ptre bending, free of shear strain, may be 
represented exactly. The complete displacement fields for this element arc given by; 


u ■ i; H.uj + (I - a*) u> -h (1 - t*) u« 
1-1 


(A. 6-6) 


T. H v. ♦ (1 - s*) V, + (I - t*) 
l-l * 


V ■ 


(A.6-7) 



The dieplacMmt wodM alileh are rapr a m tad bp t « 5 Mi i te IqaacloM 
A. 6~6t aad A. 6-7 art aaaaciacai with laMnal M gta ai of fra a dM a^ art 
illuatratad la MpiM A. 6-5. Hmm ii^laeaMBta ara Mi wa 

tef^i hy only two aoial potota aliM^ tha adip of tha a lMMt a. Siaea a 
qMiratic can Mt Miqualy iafiaai hf tmo poUita, the it^ lacMMt fiaU ' 
al^4 a coaaon ^iga of two elaaMta geoarally not equal} hanca , tiM iti^laeaBMta 
era referrai to m ineoi^tihla. 

Since me mniitim for nonotmic convarganea ia that eonpatihility anat 
axiat batwam tha alaamta. Since iiaplacaamta in gamral era not coapatibla» 
araotmic emvargmee ia not Maurei. However, ainca the incoapatibla noim 
Mtiafy tlM raquiraaenta for pure banding, m inprovad mlutim ia obtainad. 

BODHDARY COHDITIOMS 


The boundary condition at the support are illustrated in figure A. 6-3* 
Point 0 is assumed to be fixed, while other nodal points along 
lines OA and Ofi are assumed to be on roller ^^upports. -This was 
accomplished by rotating the degrees of freedoa associated with each 
of these supports into a local coordinate system that is parallel 
and normal to the roller supports. The force and dlsplaceamt 
vectoza in this coordinate systm cue: 

Cf’> • [r] (p) 

{«•> • [r] {«} 


where {P) and {6) represent the force and displacement vector 
respectively. Tlie prime indicates the local coordinate systM, 
and w is a rotation matrix as follows: 



where 9 is the angle between the glbbal x axis and the local x axis 


TIm ntrtx M.uatim U ptrfonMd oa tte •iMMit atlffMaj 
■airiett oeiy. All «l«aMits In ^ tctnc^dar perfe^.^ of th« 

■truetuTM sr« cbodMd to Aotomino if eno or aoro of Iti nodal 
points nro om of the rotated supporte. Oi« only teiw In ^ 

stiftkHMs sMtrlx that are aoctilfled are tlMoe tta^ teiiNMfond to tlw 

rotated suK>ort. 

SOLUflOM AUBOmmi 

The eon^ter proeran developed used n blo<dd.n8 al^ritto tar an out 
of core equation solver that was presented by Lestingi and Prachuktve £383 
^is approadi peisilts a very large nvad)er of equations to be solved 
on a cOTputer having liadted core capacity, ^is is beeame firing 
the asse!d>ly and decoopcsition atepa only a portion of the atrueture 
atiffheaa natrix haa to be in core. The balance can be atored on a 
disk. Similarly, during the back substitution step, only a portion 
of the decoapoaed stiffness matrix must be in core. 

It should be noted taat this solutitm algorithm takes advantage of 
the symmetzy and banded nature of the stlfftms matrix, only 
tenaa atored are the main diagonal and tte upper triangular porti<»x 
vitMn the upper band vldte. Therefore, aU of the sero teiM outside 
the band vidth are not processed. 

Since t lO is an ort.hono«aal matrix 

be stiffness S(],u&tlOD In the global coordinate system 

<p) > [k] {«> 

becomes tP') " C*] [k3 {8') 

in the local coordinate system. Ibe stiffbess mstrlx in the locsl 
astern bec^Ms 
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APPENDIX 7 
DATA INPUT 


The gear data set is input via the RAMEUST arrays defined 
in the main program. NVberical data maybe input without format statmaents* 
and fields are generated as req\iired. The input variables required, 
along with their respective NAMELIST headings «u:*e: 

/COHTRL/ 

INPUT - alj^anuBeric code used to designate type of 
input data 

•ENGL‘ - English (lbf,in. ,sec. ) 

'SI* -metric (newtons , mm, sec. ) 

OUTHJT - alphanumeric code used to designate output; 
codes used are same as for input 

IPLOT - 0; tabiilate all dynamic results 

1; plot frequency response starting at time = 0 
2; plot steady - state frequency response 

MODF - alphanumeric code used to designate whether 
or not larofile modifications are input 

'HO* - no modifications 

'YES' - modifications listed under /PRFDEF/ 

NTYPE - 1; static analysis only 

2; static and dynamic analysis only 
3; finite element analysis based on static loads 
finite element analysis based on dynamic loads 

FELGR - 1; finite element analysis based on maximm dynamic 
load experienced 

2; finite element analysis based on maximum dynamic 
load applied at pitch point 
3; finite element analysis based on maximum dynamic 
bending mcanent 

finite element analysis based on maximvm static 
bending moment 

/PHYPAR/ (two data points (one for each gear) required per variable) 

£ . Young's modulus 

PR - Poisson's ratio 

GAMA - specific weight 

* JG - polar moment of Inertia 



/GSIPAR/ 


' DP - dlaoetnl pltcli C&igllsh input 
M - gear suxlule (i^tric noduel 

DEL15» - backlaah 

TIH - input torque 

RPMm - input RPM 

ZETi^ - daaping coeffieieitt of shaft 
ZETAG - daaqiing coeffieicmt between gear teeth 
PHID - |n>e88\tre angle (degrees) 

CBD - driver bearing dan^ng coefficient 

CBl - pinion bearing dteq^lng coeffici«it 

CB2 - driven gear bearing dai^lng coefficient 
CBL - load bearing damping coefficient 

• JD - mass ncxaent of inertia of driver 

* JL - mass moaent of inertia of load 

• KDS - torsitmal S3>ring stiffness of driving shaft 

* KLS > torsional spring stiffness of load shaft 

• LDS - length of drive shaft 

* LLS - length of load shaft 


/GEOPAR/ (two data points required per variable) 


TG - number of gear teeth 

AD ' addendum 

WD - working depth 

GRRF - fillet radius of basic rack 

* Rl - hub radius 

FW - face width 

/Timu-U 

NNODE - even number of profile points used in mesh generation 
HODES - index number of those profile points used in mesh 
generation 

NGEAE - 1; stress analysis done on pinion 

2; stress analysis done on driven gear 
3; stress analysis done on gear set 

* optional, if no value entered. Program will generate values 
eis shown at the end of this section. 


The gear thooth profile can also be modified to simulate tip 
relief or undercutting. SinusoidsH errors can be introduced, as well as 
pits, to simulate involute errors due to manufacturing and surface damage, 
respectively. These modifications are introduced in /PRFDEF/ namelist. 

If MODF“ HO, /PRPDEF/ need not be included in the data card set. 

/PRPDEF/ (two data points required per variable) 

PAIW - parabolic tip modification 

ffTTO - straight line tip modification 

u.iTH - roll angle of tip modification 



I 

i 



PABM - part.bolio ^ttoa aodifloatcm 

STEM • stralc^t line bottoa aodificatlon 

RASM • roll oagle of bottoa modification 

PER - amplitude of simiTOidal error 

PAP •> phase angle of sinusoidal error 

CYC - nuBter of eyeles of sinusoidal errors 

IPIT - profile coordinate points over which pit occurs 

DEEP - depth of pit 

Q - radius to top of undercut. Fig. A. 3*4 autoauitlcally 

calculated, if not given. 

Use of the HAMEUST arrays offers a simple, unformatted 
means of Inputtlzxg data and is convenient for looping more than one data 
set. After the initial data set, subsequent data sets need Just to 
input revl8ic»is. If a later HAMELIST array contains no revisions, only 
a card with the array heading and ending need be submitted. Unlisted 
variables default to the previous values. Examples of input data card 
sets illustrate the following NAMELIST data card format (Figure A. 7-1); 

1. Column one is blank. 

2. 'A* is used to signify new NAMELIST array. 

3. 'A' is followed by the NAMELIST name. 

U. A blank separates the NAMELIST name and the first vsoriable 
name. Subsequent variables are separated by conanas. 

3. niere are two methods for defining the two elasent variables. 

The elements are defined in the order they are to be entered 
in the variable and separated by coomais, i.e., TG=32,96 defines 
TG(1)=32, and TG(2)®96. If both elements are equal, they may 
be entered by listing the number of identical values, the 
multiplication symbol, and then the value itself, i.e., AD=2*0.125, 
defines AD(l)*0.125 and AD(2)«0.125. 

6. The leist listed array value is followed by a blank and then the 

symbol from column 2 is repeated. The word END imnediately follows 
the symbol and signifies the end of that array. 

The program has the capability for accepting either SI or 
English gecu: input data and has options to print the results in either 
SI or English units. Input and output do not necessarily have to be 
of the same regime, i.e., SI output can be obtained from English input and 
vica-versa. Data sulmtitted under the *0JGL* code should be in pounds- 


•7 


forc«t ineh«s* aiul seconds. The data sutsdtted un^ar the 'SI* code 

should be In newtons, ml lllaeters, and seconds* The only oceeption to 

3 

this Is t}^ density vali^ under the 'SI' code should be in kg/n . 


PBOCffiAM DEFAULT VALUES 


Rl(l) - ( 8.0 • TIH)/ICPI • TAIMAX) ** ( 1 / 3)3 in- 
RI(2) » (8.0 * T0UT)/lPI • TADMAX) •• (1/3)1 in. 

TAWttX ■ 10,000 pBi 

JG(1) - 0.5 * GAMA(l) "PI * PH(1) » (RPCl •• U)/386. in.-lbf -sj/radlan 
JG(2) * 0.5 • Cjm{2) • PI • FM(2) * (RPC2 *• i*)/386. In.-lbf -sVradlan 

JD * 0.5 * JG(l) in. -Ibf-S^/radian 

JL » 0.5 • JG(2) in- Ibf-S /radian 

LDS » 6 in. 

LLS « 6 in. 

KDS = PI • (2.* RI(1) •* U) * [E/(2.* (l + PR)J)/(32.» LDS) in.-lbf /radian 
KLS » PI • (2.» RT.(2) ** U) • {E/(2.* (l + PR))]/(32-* US) in.-lbf /radian 


The listing of the program could be obtained by contacting the Project 
Manager at NASA Lewis Research Center. 


^C^MTC:<. 1HP'JT»*C'«GL* . •TJTPLIT»*GNGL' riPL0Tia,rtn])f*'Na*VNfyPE-4.rELGR-l S.END 
G*3«30C's.PP»'i;*f*.J55,'ja(1<^»2«0.2SO» JG»?«0'»<tHD 
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3,cB2=1.5.CBi.^1.5 &GMP 
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Figure A. 7-1 SAMPLE DATA CARD SET (W/0 PROFILE MODIFICATIONS OR ERRORS) 
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Figure A. 7-1 SAMPLE DATA CARD SET (WITH PROFILE ERRORS) ORIGINAL PAGE IS 
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SAMPLE OUTPUTS 


In the static and dynamic portions the following basic groups of 
parameters are calculated: 

a. Gear tooth and mesh deflections and the accompanying gear mesh and 
gear tooth pair stiffness. 

b. Load distribution among the contacting tooth pairs. 

c. The sliding velocity, the maximum Hertz contact pressure, and the 
sliding velocity-hertz pressure product (PV) along the tooth profiles. 


These parameters can be printed out In tabular form or plotter as Individual 
graphs on a line printer for both static and dynamic conditions. 


Some of theses plots are Illustrated for 32-96 T, 20° PA, 8 DP, standard 
full depth tooth form, gear pair, HSFsy.S, transmitting a torque equivalent 
to a normal load of 175 N/mm (1000 Ib/ln) at 2000 rpm. CR,= 1.758, CR » 
2.123. Gears are without profile modifications and errors. L = .05, c = .005 


Figure A. 8-1 
Figure A. 8-2 
Figure A. 8-3 
Figure A. 8-4 
Figure A. 8-5 
Figure A. 8-6 


Gear Tooth Deflections 

Gear Mesh and Gear Tooth Pair Stiffness 

Load Distribution 

Hertz Contact Stress 

Sliding Velocity 

PV: Hertz Contact Stress - Sliding Velocity Product 


Firues A. 8-1 through A. 8-6 are for the pinion (Gear 1) for equivalent static 
load condition at 2000 rpm. Similar plots can be generated for Gear 2, and 
as well as for the dynamic load conditions. The program accomodates the 
AGMA and metric gears, as well as the English and SI units. 


The dynamic simulation results can be also given In both the tabular and 
graph forms. For example. 


Figure A. 8-7 
Figure A. 8-8 
Figure A. 8-9 
rioiire A. 8-10 


Gear Mesh Dynamic Load 

Tabulation of Dynamic Load Factors 

Dynamic Load Between Contacting Gear Tooth Pair 

Dynamic Hertz Stress 
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Figure A. 8-1 GEAR TOOTH DEFLECTIONS 
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Figure A. 8>2 GEAR MESH AND GEAR TOOTH PAIR STIFFNESS 
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Figure A. 8-4 HERTZ CONTACT STRESS 
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Figure A. 8>6 PV: HERTZ COHTACT STRESS 
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^ieure A. 8-7 GEAR MESH DYNAMIC LOAD 
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THE INFORMATION IN THIS ABLE COMES FROM AN ANALYSIS OF THE NUMERICAL INTEORATION 
SEQUENCE COVERING THE LAS 7 PASSAGE OF A TOOTH PAIR THROUGH THE CONTACT ZONE. 
POSITION 1 CORRESPONDS TO THE STARTING POINT OF CONTACT WHILE POSITION TOO 
CORRESPONDS TO THE END POINT OF CONTACT. THE FOLLOWING SYMBOLS ARE UTILIZED IN 
THIS TABLE: 
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Figure A. 8-9 DYNAMIC LOAD BETWEEN CONTACTING GEAR TOOTH PAIR 



